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Abstract 

Distributed abstract programs are a novel class of distributed optimization problems where (i) the number of variables is much 
smaller than the number of constraints and (ii) each constraint is associated to a network node. Abstract optimization programs 
are a generalization of linear programs that captures numerous geometric optimization problems. We propose novel constraints 
consensus algorithms for distributed abstract programs: as each node iteratively identifies locally active constraints and exchanges 
them with its neighbors, the network computes the active constraints determining the global optimum. The proposed algorithms 
are appropriate for networks with weak time-dependent connectivity requirements and tight memory constraints. We show how 
suitable target localization and formation control problems can be tackled via constraints consensus. 

Index Terms 

Distributed optimization, linear programming, consensus algorithms, target localization, formation control. 

I. Introduction 

This paper focuses on a class of distributed optimization problems and its application to target localization and formation 
control. Distributed optimization and computation have recently received widespread attention in the context of distributed 
estimation in sensor networks, distributed control of actuator networks and consensus algorithms. An early established reference 
on distributed optimization is J4), whereas a non-exhaustive set of recent references includes |5|, |6|, [7], |8|. We consider a 
distributed version of abstract optimization problems. Abstract optimization problems, sometimes referred to as abstract linear 
programs or as LP- type programs, generalize linear programming and model a variety of machine learning and geometric 
optimization problems. Examples of geometric optimization problems include the smallest enclosing ball, the smallest enclosing 
stripe and the smallest enclosing annulus problems. Early references to abstract optimization problems include (9), ifTOl . ifTTI . 
In this paper we are interested in abstract optimization problems where the number of constraints n is much greater than 
the number of constraints S that identify the optimum solution (and where, therefore, there is a large number of redundant 
constraints). For example, we are interested in linear programs where n is much greater than the number of variables d 
(in linear programs, S — d). Under this dimensionality assumption, we consider distributed versions of abstract optimization 
programs, where n is also the number of network nodes and where each constraint is associated to a node. We consider 
processor networks described by arbitrary, possibly time-dependent communication topologies and by computing nodes with 
tight memory constraints. After presenting and analyzing constraints consensus algorithms for distributed abstract optimization, 
we apply them to target localization in sensor networks and to formation control in robotic networks. 

The relevant literature is vast; we organize it in three broad areas. First, linear programming and its generalizations, including 
abstract optimization, have received widespread attention in the literature. For linear programs in a fixed number of variables 
subject to n linear inequalities, the earliest algorithm with time complexity in 0(n) is given in [12|. An efficient randomized 
algorithm is proposed in [9], where a linear program in d variables subject to n linear inequalities is solved in expected time 
0(d 2 n + e°^ dlogd ' > ); the expectation is taken over the internal randomizations executed by the algorithm. An elegant survey 
on randomized methods in linear programming and on abstract optimization is |fl3l ; see also IfTTI . Ifl4l . The survey |fl5l , 
see also fl6l . discusses the application of abstract optimization to a number of geometric optimization problems. Regarding 
parallel computation approaches to linear programming, linear programs with n linear inequalities can be solved IfTTI by n 
parallel processors in time 0((loglog(n)) d ). However, the approach in [17|, see also references therein, is limited to parallel 
random-access machines, where a shared memory is readable and writable to all processors. Other references on distributed 
linear programming include |18|, |19|. 
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A second relevant literature area is distributed training of support vector machines (SVMs). A randomized parallel algorithm 
for SVM training is proposed in ll20l by using methods from abstract optimization and by exploiting the idea of exchanging 
only active constraints. Along these lines, ll2D extends the algorithm to parallel computing over strongly connected networks, 
ll22l contains a comprehensive discussion of SVM training via abstract optimization, and [23 1 applies similar algorithmic ideas 
to wireless sensor networks. The algorithms in [20|, [21 1, independently developed at the same time of our works JTJ, |2), 
[3|, differ from our constraint consensus algorithm in the following ways. First, the number of constraints stored at the nodes 
grows at each iteration so that both the memory and the local computation time at each node may be of order n. Second, our 
algorithm is proposed for general abstract optimization problems and thus may be applied to a variety of application domains. 
Third, our algorithm exploits a novel re-examination idea, is shown to be correct for time-varying (jointly strongly connected) 
digraphs, and features a distributed halting condition. 

As third and final set of relevant references, we include a brief synopsis of recent progress in target localization in sensor 
networks and formation control in robotic networks. The problem of target localization has been widely investigated and recent 
interest has focused on sensors and wireless networks; e.g., see the recent text [24|. In this paper we take a deterministic 
worst-case approach to localization, adopting the set membership estimation technique proposed in 11251 . A related sensor 
selection problem for target tracking is studied in 11261 . Regarding the literature on formation control for robotic networks, 
an early reference on distributed algorithms and geometric patterns is [27|. Regarding the rendezvous problem, that is, the 
problem of gathering the robots at a common location, an early reference is ll28l . The "circle formation control" problem, i.e., 
the problem of steering the robots to a circle formation, is discussed in ll29l . The references [ 30 1, [31], 11321 are based on, 
respectively, control-Lyapunov functions, input-to-state stability and cyclic pursuit. 

The contributions of this paper are twofold. First, we identify and study distributed abstract programming as a novel class 
of distributed optimization problems that are tractable and widely applicable. We propose a novel algorithmic methodology, 
termed constraints consensus, to solve these problems in networks with various connectivity and memory constraints: as each 
node iteratively identifies locally active constraints and exchanges them with its neighbors, the globally active constraints 
determining the global optimum are collectively identified. A constraint re-examination idea is the distinctive detail of our 
algorithmic design. We propose three algorithms, a nominal one and two variations, to solve abstract programs depending on 
topology, memory and computation capabilities of the processor network. We formally establish various algorithm properties, 
including monotonicity, finite-time convergence to consensus, and convergence to the possibly-unique correct solution of the 
abstract program. Moreover, we provide a distributed halting conditions for the nominal algorithm. We provide a conservative 
upper bound on the completion time of the nominal algorithm and conjecture that the completion time depends linearly on 
n (i.e., the number of constraints and the network dimension). Next, we evaluate the algorithm performance via a Monte 
Carlo probability-estimation analysis and we substantiate our conjecture on stochastically-generated sample problems. Sample 
problems are randomly generated by considering two classes of linear programs, taken from l33l . and three types of graphs 
(line-graph, Erdos-Renyi random graph and random geometric graph). 

As a second set of contributions, we illustrate how distributed abstract programs are relevant for distributed target localization 
in sensor networks and for formation control problems, such as the rendezvous problem and the line or circle formation 
problems. Specifically, for the target localization problem, we design a distributed algorithm to estimate a convex polytope, 
specifically an axis-aligned bounding box, containing the moving target. Our proposed eight half-planes consensus algorithm 
combines (i) distributed linear programs to estimate the convex polytope at a given instant and (ii) a set-membership recursion, 
consisting of prediction and update steps, to dynamically track the region. We discuss correctness and memory complexity 
of the distributed estimation algorithm. Next, regarding formation control problems, we design a joint communication and 
motion coordination scheme for a robotic networks model involving range-based communication. We consider formations 
characterized by the geometric shapes of a point, a line, or a circle. We solve these formation control problems in a time- 
efficient distributed manner combining two algorithmic ideas: (i) the robots implement a constraints consensus algorithm 
to compute a common shape reachable in minimum-time, and (ii) the network connectivity is maintained by means of an 
appropriate standard connectivity-maintenance strategy. In the limit of vanishing robot displacement per communication round, 
our proposed move-to-consensus-shape strategy solves the optimal formation control tasks. 

Paper organization: The paper is organized as follows. Section [TT] introduces abstract optimization problems. Section [TIT] 
introduces network models. Section IV contains the definition of distributed abstract program and the constraints consensus 
algorithms. Section [V] contains the Monte Carlo analysis of the time-complexity of the constraints consensus algorithm. 
S ections [VT| and \VU\ contain the application of the proposed constraints consensus algorithms to target localization and formation 
control. 

Notation: We let Z>o, Z>o, and M > o denote the natural numbers, the non-negative integer numbers, and the positive real 
numbers, respectively. For r 6 R>o and p £ M d , we let B(p,r) denote the closed ball centered at p with radius r, that is, 
B(p,r) = {q e R d | \\p - q\\ 2 < r}. For d e Z >0 , we Od and Id denote the vectors in M. d whose entries are all and 1, 
respectively. Similarly, we let oo^ and — oo^ the vectors with d entries 00 and —00, respectively. For a finite set A, we let 
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card(A) denote its cardinality. For two functions f,g : Z>o — * M>o, we write f(n) £ O(g) (respectively, f(n) £ £l(g)) if 
there exist N £ Z>o and c £ M > o such that f(n) < cg(n) for all n > N (respectively, f(n) > cg(n) for all n > N). Finally, 
we introduce the convention that sets are allowed to contain multiple copies of the same element. 

Given S C M. d and p £ R d , let dist(p, S) denote the distance from p to S, that is, dist(p, S) — inf se s \\p — s\\. For distinct 
pi £ R d and p2 £ Mr, let £(pi,p2) be the line through pi and pi- In what follows, a set of distinct points {pi, . . . ,p n } C R d , 
n > 3, is in stripe-generic position if, given any two ordered subsets (p a ,Pb,Pc) and (p a ,Pp,P-y), either dist(p a , £(pb,p c )) ^ 
dist(p a ,£(pp,p 7 )) or (p a ,p b ,p c ) = (p a ,pp,p 7 ). 

II. Abstract optimization 

In this section we present an abstract framework lfT31 . |[T4l that captures a wide class of optimization problems including 
linear programming and various machine learning and geometric optimization problems. Abstract optimization problems are 
also known as abstract linear programs, generalized linear programs or LP-type problems. 

A. Problem setup and examples 

We consider optimization problems specified by a pair (H,<f>), where H is a finite set, and : 2 H — > $ is a function^] with 
values in a linearly ordered set ($, <); we assume that $ has a minimum value — oo. The elements of H are called constraints, 
and for G C H, 0(G) is called the value of G. Intuitively, 0(G) is the smallest value attainable by a certain objective function 
while satisfying the constraints of G. An optimization problem of this sort is called an abstract optimization program if the 
following two axioms are satisfied: 

(i) Monotonicity: if F C G C H, then 0(F) < cf)(G); 

(ii) Locality: if F C G C H with -oo < 0(F) = cf)(G), then, for all h £ H , 

0(G) <0(Gu{/i}) cb(F) <0(FU{/i}). 

A set £> C i? is minimal if 0(-B) > 4>{B') for all proper subsets B' of -B. A minimal set B with — oo < 4>(B) is a basis. 
Given G C if, a foasis of G is a minimal subset B C G, such that — oo < 0(B) = 4>(G). A constraint /i is said to be violated 
by G, if 0(G) < 0(GU{h}). 

A solution of an abstract optimization program (if, 0) is a minimal set C if with the property that (p(Bjj) = <t>{H). 
The combinatorial dimension 5 of (7?, 0) is the maximum cardinality of any basis. Finally, an abstract program is called basis 
regular if, for any basis with card(B) = S and any constraint h £ H, every basis of B U {h} has the same cardinality of B. 
We now define two important primitive operations that are useful to solve abstract optimization problems: 

(i) Violation test: given a constraint h and a basis B, it tests whether h is violated by B; we denote this primitive by 
VLol(B,h); 

(ii) Basis computation: given a constraint h and a basis B, it computes a basis of B U {h}; we denote this primitive by 

Basis(_B, h). 

Example II.l (Abstract framework for linear programs) We recall from [ 1 3 1 how to transcribe a linear program into an 
abstract optimization program. A linear program (LP) in x £ M. d is given by 

T 

mm c x 
subject to afx < bi, i £ {1, . . . , n} 

where d £ Z>o is the state dimension, c £ R d characterizes the linear cost function to minimize, and a, £ M. d and bi £ M. 
describe n £ Z>o inequality constraints. In order to transcribe the LP into an abstract program, we need to specify the 
constraint set H and the value 0(G) for each G C H, The constraint set H is simply the set of half-spaces {hi, ... , h n }, 
where hi — {x £ Mr | afx < bi}. Defining the value function in order to satisfy the monotonicity and locality axioms is 
more delicate: if $ = (R, <) and 0(G) is the minimum of c T x subject to the constraint set G, then the locality axiom no 
longer holds (see Section 4 in lfT3ll for a counterexample). A correct choice is as follows: let ($, <) be the set M. d with the 
lexicographical order^ and define 0(G) = vq, where vq £ R d is the (unique) lexicographically minimal point x minimizing 
c T x over the constraint set G, when it exists and is bounded. If the problem is infeasible (the intersection of the constraints 
in G is empty), then vq = oo^. If the problem is unbounded (no lexicographically minimal point exists), then vq — — oo^. If 
vq is finite, then a basis of G is a minimal subset of constraints B C G such that vb = vq. It is known [9| that the abstract 
optimization program transcription of a feasible LP is basis regular and has combinatorial dimension d. A constraint h £ H 
is violated by G if and only if vq < «gu{/i}- ^ 

'Given a set H, the set 2 H is the set of all subsets of H 

2 In the lexicographic order on K 2 , we have (x\, yi) < (x2, 1/2) if and only if %i < ^2 or (xi = X2 and y\ < yi)- 
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Example II.2 (Abstract optimization problems in geometric optimization) We present three useful geometric examples, 
illustrated in Figure [T] 

(i) Smallest enclosing ball: Given n distinct points in M. d , compute the center and radius of the ball of smallest volume 
containing all the points. This problem is J9) an abstract optimization program with combinatorial dimension d+1. 

(ii) Smallest enclosing stripe: Given n distinct points in M 2 in stripe-generic positions, compute the center and the width 
of the stripe of smallest width containing all the points. In the Appendix we prove (for the first time at the best of our 
knowledge) that this problem is an abstract optimization program with combinatorial dimension 5. 

(iii) Smallest enclosing annulus: Given n distinct points in M. 2 , compute the center and the two radiuses of the annulus 
of smallest area containing all the points. This problem is [9| an abstract optimization program with combinatorial 
dimension 4. 




Fig. 1. Smallest enclosing ball, stripe and annulus 

In all three examples, the violation test and the basis computation primitives amount to low-dimensional geometric problems 
and are more or less straightforward. Numerous other geometric optimization problems can be cast as abstract optimization 
programs as discussed in 0, QU), E3, E] □ 

We end this section with a useful lemma, that is an immediate consequence of locality, and a useful rare property. 

Lemma II.3 For any F and G subsets of H, <j)(FUG) > 4>{F) if and only if there exists g £ G such that cj)(FU{g}) > <j>{F). 

Proof: If there exists g £ G such that <j>(F U {g}) > 4>(F), then by monotonicity 4>(F UG)> <f>{F U {g}) > <j>{F). For 
the other implication, assume G = {g±, . . . , gk] for some k £ Z>o, and define Gi := {gi, . . . , gi} for i £ {1, . . . , k}. We may 
rewrite the assumption <j>(F U G) > <f>(F) as <p(F U G fc _i U {gk}) > <}>(F). If 4>{F U G fc _i) = 4>{F), then the locality axiom 
implies <p(FU {gk}) > 4>(F) and the thesis follows with g = gk- Otherwise, the same argument may be applied to Gk-i- The 
recursion stops either when (j>(FU Gi) — <fi(F) (and the thesis follows with g = gi+i) for some i or when <f>(FU G\) > 4>(F) 
(and the thesis follows with g = gi). ■ 
Next, given an abstract optimization program (H,uj), let Bq denote the basis of any G C H. An element h of Bh is 
persistent if h £ Bq for all G C H containing h. An abstract optimization program (H, u>) is persistent if all elements of Bh 
are persistent. The persistence property is useful, as we state in the following result. 

Lemma II.4 Any persistent abstract optimization program (H, u>) can be solved in a number of time steps equal to the 
dimension of H. 

Proof: Let H = {hi, . . . , h n }. Set B — {hi,..., fig} and then update B = Bas±s(B,hk) for k = 5 + 1, . . . , n. 
Because of persistency, each h £ Bh is added to B once it is selected as hk and is not removed from B in subsequent basis 
computations. ■ 
Unfortunately, the persistence property is rare. Indeed, in Figure[2]we show an LP problem where the persistency property does 
not hold. In fact, it can be easily noticed that {hi, hi} is a basis for {hi, hi, /13, /14}, but {/13, /14} is a basis for {/12, ^3, h^}. 
In other words {/12} is not violated by {/13, /14}. The lack of persistency complicates the design of centralized and distributed 
solvers for abstract optimization problems. For example, in network settings, flooding algorithms are not sufficient. 

B. Randomized sub-exponential algorithm 

In the following sections we will assume that each node in the network possesses a routine capable of solving small- 
dimensional abstract optimization programs. For completeness' sake, this section reviews the randomized algorithm proposed in 
[9 |. This algorithm has expected running time with linear dependence on the number of constraints, whenever the combinatorial 
dimension 5 is fixed, and with sub-exponential dependence on the S; these bounds are proven in |9 | for linear programs and 
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Fig. 2. An LP problem that is not persistent. 

in lTP3l for general abstract optimization programs. The algorithm, called Subex_LP, has a recursive structure and is based 
on the violation test and the basis computation primitives. Given a set of constraints G and a candidate basis C C G, the 
algorithm is stated as follows: 

function Subex_LP(G, C) 
l: if G = C, then return C 
2: else 

3: choose a random h S G\C and compute B := Subex_LP(G \ {h}, C) 

4: if VLol(B, h) (that is, h is violated by B), then 

5: compute £?' := basis for B U {/i} 

6: return Subex_LP(G, Basis(B, /i)) 

7: else return B 

8: endif 

9: endif 

For the abstract optimization program (H, <fr), the routine is invoked with Subex_LP(i/, B), given any initial candidate basis 
B. 

III. Network models 

Following l34ll . we define a synchronous network system as a "collection of computing elements located at nodes of a 
directed network graph." We refer to computing elements are processors. 

A. Digraphs and connectivity 

We let Q = ({1, . . . , n), E) denote a directed graph (or digraph), where {1, . . . , n} is the set of nodes and E C {1, . . . , n} 2 
is the set of edges. For each node i of Q, the number of edges going out from (resp. coming into) node i is called out-degree 
(resp. in-degree). A digraph is strongly connected if, for every pair of nodes € {1, . . . , n} X {1, . . . , n}, there exists a 
path of directed edges that goes from i to j. A digraph is weakly connected if replacing all its directed edges with undirected 
edges results in a connected undirected graph. In a strongly connected digraph, the minimum number of edges between node 
i and j is called the distance from i to j and is denoted dist(i, j). The maximum dist(i,j) taken over all pairs is the 
diameter and is denoted diam(Cf). Finally, we consider time-dependent digraphs of the form t i— > Q{t) = ({1, . . . ,n},E(t)). 
The time-dependent digraph Q is jointly strongly connected if, for every t e Z> , the digraph \JJ= t Q(r) is strongly connected. 

In a time-dependent digraph, the set of outgoing (incoming) neighbors of node i at time t are the set of nodes to (from) 
which there are edges from (to) i at time t. They are denoted by Af out (i,t) and Afm(i,t), respectively. 

B. Synchronous networks and distributed algorithms 

A synchronous network is a time-dependent digraph Q = ({1, . . . , n}, E cmm ), where {1, . . . , n} is the set of identifiers of 
the processors, and the time-dependent set of edges E cmm describes communication among processors as follows: (i, j) is in 
E cmm (t) if and only if processor i can communicate to processor j at time t G Z>o. 

For a synchronous network Q with processors {1, . . . , n}, a distributed algorithm consists of (1) the set W, called the set 
of processor states to™, for all i € {1, . . . , n}; (2) the set A, called the message alphabet, including the null symbol; (3) 
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the map msg : W x {1, . . . , n} — > A, called the message-generation function; and (4) the map stf : W x A™ — > W, called 
the state-transition function. The execution of the distributed algorithm by the network begins with all processors in their 
start states. The processors repeatedly perform the following two actions. First, the zth processor sends to each of its outgoing 
neighbors in the communication graph a message (possibly the null message) computed by applying the message-generation 
function to the current value of vM. After a negligible period of time, the ith processor computes the new value of its processor 
state iv® by applying the state-transition function to the current value of ioM, and to the incoming messages (present in each 
communication edge). The combination of the two actions is called a communication round or simply a round. 

In this execution scheme we have assumed that each processor executes all the calculations in one round. If it is not possible 
to upper bound the execution-time of the algorithm, then one may consider a slightly different network model that allows 
the state-transition function to be executed across multiple rounds. When this happens, the message is generated by using the 
processor state at the previous round. 

The last aspect to consider is the algorithm halting, that is a situation such that the network (and therefore each processor) 
is in a idle mode. Such status can be used to indicate the achievement of a prescribed task. Formally we say that a distributed 
algorithm is in halting status if the processor state is a fixed point for the state-transition function (that becomes a self-loop) 
and no message (or equivalently the null message) is generated at each node. 

IV. Distributed abstract optimization 

In this section we define distributed abstract programs, propose novel distributed algorithms for their solutions and analyze 
their correctness. 

A. Problem statement 

Informally, a distributed abstract program consists of three main elements: a network, an abstract optimization program and 
a mechanism to distribute the constraints of the abstract program among the nodes of the network. 

Definition IV.l A distributed abstract program is a tuple (Q, (H,(f>),B) consisting of 

(i) Q = ({1, • • • ,n}, E cmm ), a synchronous network; 

(ii) (H, (p), an abstract program; and 

(iii) B : H — > {1, . . . ,n}, a surjective map called constraint distribution map that associates to each constraint one network 
node. 

If the map B is a bijection, we denote the distributed abstract program with the pair {Q, (H, (p)). A solution of (G, {H, <ft), B) 
is attained when all network processors have computed a solution to {H, <f>). 

Remark IV.2 The most natural choice of constraint distribution map B is a bijection; in this case, (i) the network dimension is 
equal to the dimension of the abstract optimization program and (ii) precisely one constraint is assigned to each network node. 
More complex distribution maps are interesting depending on the computation power and memory of the network processors. 
In what follows, we typically assume B to be a bijection. □ 

B. Constraints consensus algorithms 

Here we propose three novel distributed algorithms that solve distributed abstract programs. First, we describe a distributed 
algorithm that is well-suited for time-dependent networks whose nodes have bounded computation time, memory and in- 
degree. Equivalently, the algorithm is applicable to networks with arbitrary in-degree, but also arbitrary computation time and 
memory. Then we describe two variations that deal with arbitrary in-degree versus short computation time and small memory. 
The second version of the algorithm is well-suited for time-dependent networks that have arbitrary in-degree and bounded 
computation time, but also arbitrary memory (in the sense that the number of stored messages may depend on the number 
of nodes of the network). The third algorithm considers the case of time-independent networks with arbitrary in-degree and 
bounded computation time and memory. 

In all algorithms we consider a synchronous network Q and an abstract program (H, <j>) with H = {hi,- - • ,h n } and with 
combinatorial dimension 8. We define a distributed abstract program by assuming that constraints and nodes are in a one-to- 
one relationship, and we let hi be the constraint associated with network node i. Here is an informal description of our first 
algorithm. 

Constraints Consensus Algorithm: Beside having access to the constraint hi, the ith processor state contains a 
candidate basis consisting of 8 elements of H. The processor state is initialized to 8 copies of hi. At 
each communication round, the processor performs the following tasks: (i) it transmits £?M to its out-neighbors and 
acquires from its in-neighbors their candidate bases; (ii) it solves an abstract optimization program with constraint set 
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given by the union of: its constraint hi, its candidate basis flH and its in-neighbors' candidate bases; (hi) it updates 

£?W to be the solution of the abstract program computed at step (ii). 
For completeness' sake, the following table presents the algorithm in a way that is compatible with the model given in 
Section III-B The Subex_LP algorithm is adopted as local solver for abstract optimization programs. 



Constraints Consensus 

A = H 5 U {null} 
SM C H with card(SW) = 5 
SW : {//, //,! 

function tasg(B^',j) 
l: return B® 

function stf(B®,y) 

% executed by node i,,with,v,, := msg(S^.z) = B^ 
1: fltap := {M U BW U ( U i6Ml(i) 

2: return Subex_LP(iJ tmp , B®) 



Problem data: 
Algorithm: 
Message alphabet : 
Processor state: 
Initialization : 



Remark IV.3 (Constraint re-examination due to lack of persistency) In order for the algorithm to compute a correct 
solution, it is necessary that each node continuously re-examine its associated constraint throughout algorithm execution. In 
other words, step 1: of the state-transition function stf in the algorithm may not be replaced by H tmp := B<- 1 ' U ( U^-g^fi) Vj)- 



This continuous re-examination is required because of the lack of the persistency property discussed after Lemma II. 4 □ 



In the second scenario we consider a time-dependent network with no bounds on the in-degree of the nodes and on the memory 
size. In this setting the execution of the Subex_LP may exceed the computation time allocated between communication rounds. 



To deal with this problem, we introduce an "asynchronous" version of the network model described in Section III we allow a 
processor to execute message-transmission and state-transition functions at instants that are not necessarily synchronized. Here 
is an informal description of the algorithm. 

Multi-round constraints consensus algorithm Each processor has the same message alphabet, processor state, and 
initialization settings as in the previous constraints consensus algorithm. The processor performs two tasks in parallel. 
Task #1: at each communication round, the processor transmits to its out-neighbors its candidate basis B^ and 
acquires from its in-neighbors their candidate bases. Task #2: independently of communication rounds, the processor 
repeatedly solves an abstract optimization program with constraint set given by the union of: its constraint hi, its 
candidate basis Z?M and its in-neighbors' candidate bases; the solution of this abstract program becomes the new 
candidate basis fiW. The abstract program solver is invoked with the most-recently available in-neighbors' candidate 
bases and, throughout its execution, this information does not change. 
In the third scenario we consider a time-independent network with no bounds on the in-degree of the nodes. We suppose that 
each processor has limited memory capacity, so that it can store at most D constraints in H . The memory is dimensioned so 
as to guarantee that the abstract optimization program is always solvable during two communication rounds (e.g., by adopting 
the Subex_LP solver). The memory constraint is dealt with by processing only part of the incoming messages at each round, 
and by cycling among incoming messages in such a way as to process all the messages in multiple rounds. 

Cycling constraints consensus algorithm The processor state contains and initializes a candidate basis B<- 1 ' as in the 
basic constraints consensus algorithm. Additionally, the processor state includes a counter variable that keeps track 
of communication rounds. At each communication round, the processor performs the following tasks: (i) it transmits 
B^ to its out-neighbors and receives from its in-neighbors their candidate bases; (ii) among the incoming messages, 
it chooses to store D messages according to a scheduled protocol and the counter variable; (iii) it solves an abstract 
optimization program with constraint set given by the union of: its constraint hi, its candidate basis B™ and the D 
candidate bases from its in-neighbors; (iv) it updates flW to be the solution of the abstract program computed at step 
(iii). 



C. Algorithm analysis 

We are now ready to analyze the algorithms. In what follows, we discuss correctness, halting conditions, memory complexity 
and time complexity. 
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Theorem IV.4 (Correctness of the constraints consensus algorithm) Let (Q, (H, <ft)) be a distributed abstract program with 
nodes and constraints in one-to-one relationship. Assume the time-dependent network Q is jointly strongly connected. Consider 
a constraint consensus algorithm in which a network node initializes its candidate basis to a constraint set with finite value. 
The following statements hold: 

(i) along the evolution of the constraints consensus algorithm, the basis value t \— > <f>(B^ (t)) at each node i G {1, . . . , n} 
is monotonically non-decreasing and converges to a constant finite value in finite time; 

(ii) the constraints consensus algorithm solves the distributed abstract program {Q, (H, if))), that is, infinite time the candidate 
basis at each node i is a solution of (H,cf>); and 

(iii) if the distributed abstract program has a unique minimal basis B^, then the final candidate basis B^' at each node i 
is equal to Bh- 

Proof: From the monotonicity axiom of abstract optimization and from the finiteness of H, it follows that each sequence 
4>(B^(t)), t £ Z>o, is monotone non-decreasing, upper bounded and can assume only a finite number of values. Additionally, 
in finite time, each node has a candidate basis that has finite value because, by assumption, one node starts with a candidate 
basis with finite value and the digraph is jointly strongly connected. Therefore, the constraints consensus algorithm at every 



node converges to a constant candidate basis with finite value in a finite number of steps. This concludes the proof of fact (i) 
In what follows, let B^\ . . . , B>- n ' denote the limiting candidate bases at each node in the graph. 



We prove fact (ii) in three steps. First, we proceed by contradiction to prove that all the nodes converge to the same value (but 
not necessarily the same basis). The following fact is known: if a time-dependent digraph is jointly strongly connected, then 
the digraph contains a time-dependent directed path from any node to any other node beginning at any time, that is, for each 
t £ Z>o and each pair there exists a sequence of nodes tx,.-.,lk and a sequence of time instants t\,... , ifc+i £ Z>o 

with t < t\ < ■■ ■ < ffc+i, such that the directed edges {(i,£i), (^1,^2), ■ ■ ■ , (^k,j)} belong to the digraph at time instants 
{ii, . . . , tk+i}, respectively. The proof by contradiction of a closely related fact is given in [35 Theorem 9.3]. Now, suppose 
that at time to all the nodes have converged to their limit bases and that there exist at least two nodes, say i and j, such that 
<f>(B®) ^ cj)(B^). For t = t + l, for every k\ £ A/" ut(£> *o + 1), no constraint in B^ violates B^ kl \ otherwise node k\ would 
compute a new distinct basis with strictly larger value, thus violating the assumption that all nodes have converged. Therefore, 
0(£?M) < 0(ij[ fe i]). Using the same argument at t = to + 2, for every ki £ A/" ou t(fci, to + 2), no constraint in B^- kl ^ violates 
B^- k2 \ Therefore, <f)(B^) < tfi(B^ kl ^) < cj)(B^). Iterating this argument, we can show that for every S > 0, every node k, 
that is reachable from i in the time-dependent digraph with a time-dependent directed path of length at most S, has a basis 
£?[ fe l such that (f>{B^) < (j)(B^). However, because the digraph is jointly strongly connected, we know that there exists a 
time-dependent directed path from node i to node j beginning at time to> mus showing that (f>(B^) < 4>{B^). Repeating the 
same argument by starting from node j we obtain that cj>(B^) < cf>(B^). In summary, we showed that <f){B^) — (f>{B^'), 
thus giving the contradiction. Note that this argument also proves that, if (i,j) is an edge of the digraph U^t^(r), then no 
constraint in i violates and, therefore, cj>(B^ U B®) = <j>(B^). Also, the equality B® U B® = B® U B® implies that 
there exists cf> £ R such that = <f>(B®) = <j>(B^) = <j>(B® U B®) = <j>(B^ U B®) for all i G {1, . . . , n} and (i,j) edges 
of U+~S(t). 

Second, we claim that the value of the basis at each node is equal to the value of the union of all the bases. In other words, 
we claim that 

4> = 4>(B [1] \J---\JB [n] ). (1) 

We prove equation ([TJ by induction. First, we note that <j> = <fi(B^ U B^) for any nodes i and j such that either (i,j) or 
(j,i) is a directed edge in \J±^Q{r). Without loss of generality, let us assume i = 1 and j — 2. Now assume that 

<f)(B [1] U •• ■ UB [k] ) = 4>, (2) 

for an arbitrary /c-dimensional weakly-connected subgraph Gk of U^^G(t) and we prove such a statement for a weakly- 
connected subgraph of dimension k + 1 containing Gk- By contradiction, we assume the statement is not true for k + 1. 
Assuming, without loss of generality, that node k + 1 is connected to Gk in U+™(/(t), we aim to find a contradiction with 
the statement 

(j>(B [1] U---UB [k] UB [k+1] ) > </>. (3) 
Plugging the induction assumption into equation ((3), we have 

4>{B [1] U • • • U B [k] U B [k+1] ) > 4>{B [1] U • • • U B [k] ). (4) 
with F = £?W U • • -U.BM and G = i?[ fc+1 l and noting equation Q, we conclude that there exists g £ B^ k+1 ^ 

cj)(B [1] U ■ • • U B [k] U {<?}) > <j)(B [1] U ■ ■ • U B [k] ). (5) 



From Lemma 
such that 



II. 3 
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Next, select a node p e {1, . . . , fc} such that either (p, k + 1) or (k + l,p) is a directed edge of U£^(/(t) and note that 
fiW C U • • • U and 0(i? [p ]) = <p(B^ U • • • U B^) by the induction assumption. From these two facts together with 
equation |5]), the locality property implies that 

#flMu{ fl })>0(BM). (6) 

Finally, the contradiction follows by noting: 

0(2JW) «— ^ L co„ve rg ed ^ (flW u B[fe+1]) —city ^ ^ j ^ by e q u. 10 „ Q 

This concludes the proof of equation (fT). 

Third and final, because no constraint in {hi, . . . , /i n } violates the set £?M U • • • U Si"' and because £?M U • • • U C if, 
Lemma II. 3 and equation ([T| together imply 

|Wu-UB w )=^(4 



This equality proves that in a finite number of rounds the candidate basis at each node is a solution to (if, </>) and, therefore, 



this concludes our proof of fact (ii) The proof of fact (iii) is straightforward 



Remark IV.5 (Correctness of multi-round and cycling constraints consensus) Correctness of the other two versions of 
the constraints consensus algorithm may be established along the same lines. For example, it is immediate to establish that 
the basis at each node reaches a constant value in finite time. It is easy to show that this constant value is the solution 
of the abstract optimization program for the multi-round algorithm over a time-dependent graph. For the cycling algorithm 
over a time-independent graph we note that the procedure used to process the incoming data is equivalent to considering a 
time-dependent graph whose edges change with that law. □ 



Theorem IV.6 (Halting condition) Consider a network described by a time-independent strongly-connected digraph Q im- 
plementing a constraints consensus algorithm in which a network node initializes its candidate basis to a constraint set with 
finite value. Each processor can halt the algorithm execution if the value of its basis has not changed after 2diam(CJ) + 1 
communication rounds. 

Proof: For all t G Z> and for every edge of Q, 

4>(B^(t))<4>(B^(t+l)), (7) 

because, by construction along the constraints consensus algorithm, the basis B^'(t + 1) is not violated by any constraint 
in the basis B^(t). Assume that node i satisfies B^(t) = B for all t G {to,..., to + 2 diam(CJ)}, and pick any other 
node j. Without loss of generality, set to = 0. Because of equation Q, if ki G J\f ut(i), then <f>(B^ (1)) > 0(B) and, 
recursively, if &2 € A<out(&i), then 0(f?[ fe2 l (2)) > c/)(B^ kl ^ (1)) > </>(B). Therefore, iterating this argument dist(i, j) times, 
the node j satisfies <p(B^(dist(i, j))) > 4>(B). Now, consider the out-neighbors of node j. For every k$ G A/" ou t(j), it 
must hold that (/}(B^ k ^ (dist(z, j) + 1)) > <j)(B^ (dist(i, j))). Iterating this argument dist(j, i) times, the node i satisfies 
<f>(B®(dist(iJ) +dist(j,«))) > <f>(BW(dist(i,j))). In summary, because dist(z, j) + dist (J, i) < 2diam(CJ), we know that 
<p(B^(dist(i,j) +dist(j, i))) = <f>(B) and, in turn, that 

4>(B) >#B0(dist(i,j))) > ct>(B). 

Thus, if basis i does not change for 2 diam(C?) + 1 time instants, then its value will never change afterwards because all bases 
B^\ for j G {1, . . . , n}, have cost equal to <f>(B) at least as early as time equal to diam(<?) + 1. Therefore, node i has sufficient 
information to stop the algorithm after a 2diam(Cf) + 1 duration without value improvement. ■ 
Next, let us state some simple memory complexity bounds for the three algorithms. Assume that 6 is the combinatorial 
dimension of the abstract program (H, (j)) and call a memory unit is the amount of memory required to store a constraint in 
H. Each node i of the network requires 1 + 6(1 + card(A/i„(i))) G 0(n) memory units in order to implement the constraints 
consensus algorithm and its multi-round variation, and 1 + 5(1 + D) G 0(1) in order to implement the cycling constraints 
consensus algorithm. 

We conclude this section with some incomplete results about the completion time of the constraints consensus algorithm, i.e., 
the number of communication rounds required for a solution of the distributed abstract program, and about the time complexity, 
i.e., the functional dependence of the completion time on the number of agents. First, it is straightforward to show that there 
exist distributed abstract programs of dimension n for which the time complexity can be lower bounded by f2(n). Indeed, it 
takes order n communication rounds to propagate information across a path graph of order n. On the other hand, it is also easy 
to provide a loose upper bound by noting that (i) the number of possible distinct bases for an abstract optimization program 
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with n constraints and combinatorial dimension S is upper bounded by n s , and (ii) at each communication round at least 
one node in the network increases its basis. Therefore, the worst-case time complexity is upper bounded by n s+1 . It is our 
conjecture that the average time complexity of the constraints consensus algorithm is much better than the loose analysis we 
are able to provide so far. 

Conjecture IV.7 (Linear average time complexity) Over the set of time-independent strongly-connected digraphs and dis- 
tributed abstract programs, the average time complexity of the constraints consensus algorithm belongs to 0(diam((?)). □ 

V. Monte Carlo analysis of the time complexity of constraints consensus 

This section presents a simulation-based analysis of the time complexity of the constraints consensus algorithm for stochas- 
tically generated distributed abstract programs. To define a numerical experiment, i.e., a stochastically-generated distributed 
abstract program, we need to specify (1) the communication graph, (2) the abstract optimization problem and (3) various 
parameters describing a nominal set of problems and some variations of interest. We discuss these three degrees of freedom 
in the next three subsections and perform two sets of Monte Carlo analysis in the two subsequent subsections. 

A. Communication graph models 

We consider time-independent undirected communication graphs generated according to one of the following three graph 
models. The first model is the line graph. It has bounded node degree and the largest diameter. Then we consider two random 
graphs models, namely the well-known Erdds-Renyi graph and random geometric graph. In the Erdos-Renyi graph an edge 
is set between each pair of nodes with equal probability p independently of the other edges. A known result |36| is that the 
average degree of the nodes is pn. Also, if p = (1 + e) log(n.)/n, e > 0, the graph is almost surely connected and the average 
diameter of the graph is log(n)/log(pn). Indeed, we use the probability (1 + e)log(n)/n to generate the graph. This gives 
an average node degree that is unbounded as n grows, but the growth is logarithmic and so the local computations are still 
tractable. A random geometric graph in a bounded region is generated by (i) placing nodes at locations that are drawn at 
random uniformly and independently on the region and (ii) connecting two vertices if and only if the distance between them 
is less than or equal to a threshold r > 0. We generate random geometric graphs in a unit-length square of R 2 . To obtain a 
connected graph we set the radius r to the minimum value that guarantees connectivity. 

B. Linear programs models 

In our experiments we consider stochastically-generated linear programs according to well-known models proposed in the 
literature. A detailed survey on stochastic models for linear programs and their use to study the performance of the simplex 
method is [33]. We consider standard LPs in d-dimensions with n constraints of the form 

T 

mm c x 
subject to Ax < b, 

where A G W lXd , b G K™ and c € R d are generated according to the following stochastic models. 

Model A. In this model the elements A^ and Cj are independently drawn from the standard Gaussian distribution. The 
vector b is defined by 6j = (X^=i Afj) 1 ^ 2 , i G {1, . . . ,n}. This corresponds to generating hyperplanes (corresponding to the 
constraints) whose normal vectors are uniformly distributed on the unit sphere and that are at unit distance from the origin. 
The LP problems generated according to this model are always feasible. This model was originally proposed by ll37l and is a 
special case of a class of models introduced by ll38l see ll33l for details^] 

Model B. In this model the vector c is obtained as c = A T c. The vector (b, c) G M. 2n is uniformly randomly generated in 
[0, l] 2n and A is a standard Gaussian random matrix independent of (b,c). The LP problems generated according to this model 
are always feasible. This LP model, with a more general stochastic model for (b,c), was proposed by Todd in [39| (where it 
is the model indicated as "Model 1" in a collection of three). 

3 Our Model A is the one indicated as Model N in 1331 . A slightly different model, called Model O, features hyperplanes at random distance from the 
origin. Numerical experiments show that, for large n/d, Model N problems require more iterations on average than Model O ones. Indeed, for fixed d and 
increasing n, all constraints in Model N problems are relevant, whereas many constraints in Model O problems are not. 
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C. Nominal problems and variations 

First, as nominal set of problems we consider a set of distributed abstract programs with the following characteristics: d = 4, 
n £ {20, . . . , 240}, the graphs are equal to the line graphs of dimension n, and the linear problems are generated from Model 
A. 

Second, as variations of the nominal set of problems, we generate LPs of dimension d 6 {2,3,4,5} with a number of 
constraints n € {20,40, . . . , 240}. For each value of d we generate a graph according to one of the three graph models and 
an LP according to one of the two LP models. For each configuration (dimension, number of constraints, graph model, and 
LP model), we generate different problems, we solve each problem with the constraints consensus algorithm, and we store 
worst-case and average completion time. 

Results for the nominal set of problems and for its variations are given in the next sections. 



D. Time complexity results for nominal problems 

For the nominal set of problems, we study the time complexity via the Student i-test and via Monte Carlo probability 
estimation. 

For each value of n, we perform a Student's i-test with the null hypothesis being that the average completion time divided 
by the graph diameter is greater than 1.5 — against the alternate hypothesis that the same ratio is less than or equal to that 
(at the 95% confidence level). (Note that the diameter is (n - 1).) The results for n € {200,220,240} are shown in Table [I] 
The tests show that we can reject the null hypothesis. In Figure [3] we show the linear dependence of the completion time with 
respect to the number of agents (and therefore with respect to the diameter) and provide the corresponding 95% confidence 
intervals. 

number of constraints average completion stan( j arc j deviation df t-value p-value 

time/diam 



240 1.21 0.36 99 -7.73 4.3 x 10~ 12 

220 1.16 0.31 99 -10.90 6.0 x 10~ 19 

200 1.27 0.36 99 -6.49 1.7 x 10~ 9 

TABLE I 

Student's t-test results for demonstrating the linear dependence of the completion time with respect to the diameter. Problem: 
Graph = line graph, LP model = Model A, d = 4, n e {200, 220, 240}, N mn = 100, null hypothesis: average completion time /n > 1.5. 




Next, we aim to upper bound the worst-case completion time. To do so we use a Monte Carlo probability estimation method, 
that we review from [40|. 

Remark V.l (Probability estimation via Monte Carlo) We aim to estimate the probability that a random variable is less 
than or equal to a given threshold. Let Q be a compact set and A be a random variable taking values in Q. Given a scalar 
threshold 7, define the probability ^(7) = Pr{ J(A) < 7}, where J : Q — > K is a given measurable performance function. We 
estimate p(j) as follows. First, we generate N € Z>o independent identically distributed random samples A^, . . . , A^ N \ 
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Second, we define the indicator function Ij n : Q — * {0, 1} by Ij 7 (A) = 1 if J(A) < 7, and otherwise. Third and final, we 
compute the empirical probability as 

1 N 

Next, we adopt the Chernoff bound in order to provide a bound on the number of random samples required for a certain 
level of accuracy on the probability estimate. For any accuracy e £ (0, 1) and confidence level 1 — 77 £ (0, 1), we know that 
IIpaKt) ~ P{l)\\ < e w i tn probability greater than 1 — 77 if 

A>-^log-. (8) 
2e z 77 

For e = 77 = 0.01, the Chernoff bound ^ is satisfied by N = 27000 samples. □ 

Adopting the same notation as in Remark [V. 1 1 here is our setup: First, the random variable of interest is a collection of n 
unit-length vectors (i.e., our random variable takes values in the compact space {x £ M 4 | ||x|| = 1}™). Second, the function 
J is the completion time of the constraints consensus algorithm in solving a nominal problem with constraints determined by 
A. Third, we want to estimate the probability that, for n £ {40, 60, 80}, the completion time is less than or equal to 4 times 
the diameter of the chain graph of dimension n. In order to achieve an accuracy 0.01 with confidence level 99%, we run 
A run = 27000 experiments for each value of 77 and we compute the maximum completion time in each case. The experiments 
show that for each n the worst-case completion time is less than 3.4 times the graph diameter. Therefore, we have established 
the following statement. 

With 99% confidence level, there is at least 99% probability that a nominal problem (d — 4, graph = line graph, LP 
model = Model A) with number of constraints n £ {40, 60, 80} is solved by the constraints consensus algorithm in 
time upper bounded by 4 (77 — 1). 

E. Time complexity results for variations of the nominal problems 

Next we perform a comparison among different graph models, LP models and LP dimensions. In order to compare the 
performance of different graphs we consider problems with: Graph £ {line graph, Erdos-Renyi graph, random geometric 
graph}, LP model = Model A, d = 4, 77 £ {20, . . . , 240}, N mn = 100. We compute the average completion time to diameter 
ratio for increasing values of 77. The results with the 95% confidence interval are shown in Figure |4] 
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Fig. 4. Comparing completion time to diameter ratio for distinct graphs and for increasing number of constraints n. Specifically: Graph S {line graph 
(circle), Erdos-Renyi graph (square), random geometric graph (diamond)}, LP model = Model A, d = 4, n £ {20, . . . , 240}, N mn = 100, performance = 
average. 

To compare the performance for different LP models, we consider problems with: Graph = line graph, LP model £ {Model 
A, Model B}, d — 4, 71 £ {20, . . . ,240}, A^,, = 100. The results with the 95% confidence interval are shown in Figure [5] 

Next, to compare the performance for different dimensions d, we consider problems with: Graph = line graph, LP model = 
Model A, d £ {2, 3, 4, 5}, 77 £ {20, . . . , 240}, A run = 100. The results with the 95% confidence interval are shown in Figure [6] 
The comparisons show that, the linear dependence of the completion time with respect to the number of constraints is not 
affected by the graph topology, the LP model and the dimension d. As regards the dimension d, as expected, for fixed 77 the 
average completion time grows with the dimension. Also, the growth appears to be linear for d > 3 (for d = 2 the algorithm 
seems to perform much better). In Figure [7] we plot the least square value of the completion time to diameter ratio over the 
number of agents 77 versus the dimension d. 
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Fig. 5. Comparing completion time to diameter ratio for distinct LP models and for increasing number of constraints n. Specifically: Graph = line graph, 
LP model e {Model A (circle), Model B (diamond)} d = 4, n e {20, . . . , 240}, N mn = 100, performance = average. 
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Fig. 6. Comparing completion time to diameter ratio for distinct problem dimensions d and for increasing number of constraints n. Specifically: Graph ; 
line graph, LP model = Model A, d S {2(square), 3(diamond),A(circle), 5(triangle)}, n S {20, . . . , 240}, N mu = 100, performance = average. 



VI. Application to target localization in sensor networks 

In this section we discuss an application of distributed abstract programming to sensor networks, namely a distributed 
solution for target localization. Essentially, we propose a distributed algorithm to approximately compute the intersection of 
time- varying convex poly topes. 

A. Motion and sensing models 

We consider a target moving on the plane with unknown but bounded velocity. The first-order dynamics is given by 

p(t + l)=p(t)+v(t), (9) 

where p(t) £ [x m ; n , £ max ] x [ymi n ,y m ax] C M 2 is the target position at time t € Z> and v(t) £ K 2 is unknown but satisfies 
IMI < Vmax for given u max . 

A set of sensors {1, . . . , n} deployed in the environment measures the target position. We assume that the measurement 
noise is unknown and such that, at each time instant, each sensor i measures a possibly-unbounded region of the plane, 
/i^(p(i)) C K 2 , containing the target. The set M(p(t)) = n^g/i. n }/iW (p(t)), called the measurement set, provides the 
best estimate of the target position based on instantaneous measures only. An example scenario is illustrated in Figure [8] In 
what follows, we make a critical assumption and a convenient one. First, we assume that each measured region h\ l \p{t)) is 
a possibly-unbounded convex polygon. Second, for simplicity of notation, we assume that each measured region h^(p(t)) is 
a half-plane, so that the measurement set M(p(t)) is a non-empty possibly-unbounded convex polygon equal with up to n 
edges. 
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Fig. 7. Comparing completion time to diameter ratio (least square value in n) for increasing problem dimensions d. Specifically: Graph = line graph, LP 
model = Model A, d £ {2, 3, 4, 5}, n 6 {20, . . . , 240}, iVnm = 100, performance = average. 



Fig. 8. The measurement set is the intersection of the 4 sensor measurements. 



B. Set-membership localization of a moving target 

Problem VI. 1 (Set-membership localization) Compute the smallest set E(t) C R 2 that contains the target position p(t) at 
t G Z>o and that is consistent with the dynamic model (|9) and the sensor measurements h^'(p(t)), i € {1, ... ,n}. □ 

We adopt the set-membership approach described in l25l j^]For t € Z>o, define the sets E(t\t— 1) and E(t\t) as the feasible 
position sets containing all the target positions at time t that are compatible with the dynamics and the available measurements 
up to time t — 1 and t, respectively. With this notation, the recursion equations are: 

£(0|0)=Af(p(0)), (10a) 
E(t\t - 1) = E(t - l\t - 1) + B(Q, t> max ), (10b) 
E(t\t) = E(t\t- l)n M(p(t)), (10c) 

where the sum set of two sets A + B is defined as {a + b | a € A and & G B}. Equation ( |10b[ ) is justified as follows: since 
the target speed satisfies ||w|| < u max , if the target is at position p at time t, then the target must be inside B(p, v max r) at 
time t + r, for any positive r. Equation ( | 1 0c] > is a direct consequence of the definition of measurement set. Equations ( |10a| i, 
( |10b[ ), and ( | 1 0c| > are referred to as initialization, time update and measurement update, respectively. The time and measurement 
updates are akin to prediction and correction steps in Kalman filtering. 

Running the recursion ( fT0| > in its exact form is often computationally intractable due to the large and increasing amount of 
data required to describe the sets E(t\t— 1) and E(t\t). To reduce the computational complexity one typically over-approximates 
these sets with bounded-complexity simple-structure sets, called approximating sets. For example, a common approximating 

4 In 1251 a moving vehicle localizes itself by measuring landmarks at known positions. This "self-localization" scenario leads to a mathematical model 
closely related to the one we consider. 
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E(0\0) = M(p(0)) 




/£(0|0> 
E(1|0) 




E(l|l) A/(p(l)) 



Fig. 9. Set-membership localization recursion (10) : initialization, time update and measurement update. 



set is the axis-aligned bounding box, i.e., the smallest rectangle aligned with the reference axes containing the set. If II denotes 
the projection from the subsets of M 2 onto the collection of approximating sets, then the recursion ( fT0] > is rewritten as 

£7(0|Q)=n(M(p(0))), (11a) 
E(t\t - 1) - U{E(t - l\t - 1) + B(0, tw)), (1 lb) 

E(t\t) = U(E(t\t - 1) n M (p(t))), (11c) 

where now the sets E(t\t — 1) and E(t\t) are only approximations of the feasible position sets. 

C. A centralized LP-based implementation: the eight half-planes algorithm 

In this section we propose a convenient choice of approximating sets for set-membership localization and we discuss the 
corresponding time and measurement updates. We begin with some preliminary notation. 

We let Hk be the set containing all possible collections of k half-planes; in other words, an element of Hk is a collection 
of k half-planes. Given an angle 9 € [0, 2ir{ and a set of half-planes H = {hi, . . . , hk} € Hk with fcj = {i e I 2 afx < 
bi, ||<Zi|| = 1, bi € K.}, define the linear program (H, 9) by 

max [cos(#) sin(0)] • x 
subject to afx < bi, i € {1, . . . , k}, (12) 

5^ X ^ X max , ymin — y — Vmnx • 



^min 



As in Example II. 1 transcribe (H,9) into an abstract optimization program (H,<j)g). Recall that (H,(j)g) is basis regular and 
has combinatorial dimension 2, so that its solution, i.e., the lexicographically-minimal minimum point of (H, 9), is always a 
set of 2 constraints, say h^ g and h^ e . In other words, the pair {^'j,'*^} G Ti-2 is computed as a function of an angle 
9 e [0, 27r[ and of a fc-tuple H = {hi, . . . , h k } € H k - 

Now, as collection of approximating sets we consider the set 7i 8 containing the collections of 8 half-planes. Note that the 
subset of elements {hi, . . . , h%} S 7is such that n^ =1 hj is bounded is in bijectiorj^] with the set of convex polygons with at 
most 8 edges. Additionally, for arbitrary k, we define the projection map IIlp : Hk — > Hs as follows: given H £ Hk, define 
n LP (iJ) S 7^8 to be the collection of half-planes h}3 e and hM e , for f? S {0, tt/2, 7T, 37r/2}. Note that our approximating set 
IIlp (if) contains H and is contained in the smallest axis-aligned bounding box containing H; additionally, note that IlLp(-ff) 
contains some possibly repeated half-planes because the same half-plane could be part of the solution for distinct values of 6. 

Our definition of IIlp has the following interpretation: assuming the target is known to satisfy all half-plane constraints in a 
set H, the reduced-complexity possibly-unbounded polygon containing the target is computed by solving four linear programs 



(H,6), 9 e {0, tt/2, 7r,37r/2}; see Fi g ure 10 



Finally, we review the approximated set-membership localization recursion ( [TT| . We assume hM'(t) — {x € E 2 | a^(t) T x < 
||aj|| = 1} is the half-plane containing the target measured by sensor i E {1, . . . , n} at time t € Z>o- The approximated 
feasible position sets, elements of Hs, are 

E(t\t - 1) = {hi(t\t - 1), ... , h 8 (t\t - 1)}, and 
E(t\t) = {hi(t\t),...,h 8 (t\t)}. 

Initialization: Equation ( | 1 1 a| > reads 

{Ai(o|o), . . . ,M0|0)} - n LP ({/j [11 (o), . . . , h^(o)}). 



5 Indeed, any convex polygon may be defined as the intersection of a finite number of half-planes; this definition is referred to as the H-representation of 
the convex polygon. 
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Time update: Assume hi(t\i) = {x G M. 2 \ ai(t) T x < bi(t)}, that is, hi(t\t) is characterized by the coefficients (di(t), 6j(i)). 
Since the target speed satisfies ||u|| < v max , at instant t + r the target is contained in the half-planes 

hi(t + r|t) = {.x e M 2 | a^tfx < bi(t) + wr}. 

Therefore, the time update consists in defining each hi(t + l|t) to be (cii(t), bi(t) + v max ); we refer to this operation as to a 
time-translation by an amount u max of the half-plane ft,j(i|t). This time-update operation is equivalent to equation (JTTbJ and 
does not explicitly require an application of the projection IIlp. 

Measurement update: Equation \\\c\ reads 

{h 1 (t\t),...,h 8 (t\t)} = U hP (H(t)), 
where the collection of constraints at time t is 

H{t) = {ftx(t|i - 1), ... , h 8 (t\t - 1)} U {fcN(t), . . . , ftW(t)}. 

Remark VI.2 (i) During each iteration of the localization recursion, the time update step requires 8 sums and the measure- 
ment update steps requires the solution of 4 linear programs in 2 variables and n + 8 constraints, 
(ii) Similar localization algorithms arise by selecting £ > 3 and by solving £ LP at each iteration parametrized by 9 E 
{0, 2n/£, . . . , 2(£ — l)w/£}. Larger values of £ lead to tighter approximating polygons. □ 



D. A distributed eight half-planes algorithm 

We consider a scenario in which the sensors measuring the target position also have computation and communication 



capabilities so that they form a synchronous network as described in Section III Let Q in = ({1, . . . ,n}, E sn ) be the undirected 
communication graph among the sensors {1, ,..,n)\ assume Q is connected. Assume the sensors communicate at each time 
t £ Z>o and perform measurements of the target at unspecified times in Z>o (communication takes place at higher rate than 
sensing). For simplicity, we assume the first measurement at each node happens at time 0. 

We aim to design a distributed algorithm for the sensor network to localize a moving target. The idea is to run local set- 
membership recursions (with time and measurement updates) at each node while exchanging constraints in order to achieve 
constraints consensus on a set-membership estimate. Distributed constraint re-examination is obtained as follows: at each time, 
each node keeps in memory the last m measurements it took and, after an appropriate time-update, re-introduces them into 
the n LP computation. We begin with an informal description. 

Eight Half-Planes Consensus Algorithm: The processor state at each processor i contains a set H^ timal of 8 candidate 
optimal constraints and a set -ff se , lsed containing the last m measurements, for some m > 0. These sets are initialized 
to the first sensor measurement. At each communication round, the processor performs the following tasks: (i) it 
transmits H^ timal to its out-neighbors and acquires from its in-neighbors their candidate constraints; (ii) it performs a 
time-update, that is, a time-translation by an amount w max , of all candidate optimal, measured and received constraints; 
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(iii) it updates the set of measured constraints if a new measurement is taken; and (iv) it updates H^ timal to be the 
projection IIlp of all candidate optimal, measured and received constraints. 
Next we give a pseudo-code description. 

Problem data : A network Q sn of sensors that measure half-plane constraints 



Algorithm: 
Message alphabet : 
Processor state: 

Initialization : 

function msg((ff j timal , H^J, j) 



Eight Half-Planes Consensus 

A = TLs U {null} 



^optimal 6 ^8 

t= T-f 

-"sensed t ri m 



for some m > 



H, 



»] 



optimal 



H. 



sensed 



:={/*W(0), 
= {/iW(0), 



.,^(0)} 
,fcM(0)} 



l: return H 



optimal 



function stf((^ timal ,ff s l l ed ),y 
% executed by^ node i, with yj := H t 

2 
3 
4 
5 



time-translate by an amount i; max an constraints constraints in -ff sensed , flHteiai' anc ' UjeA/h(i)% 
if a new measurement is taken at this time, then 

add it to -ff sensed ; drop oldest measurement from -ff sensed 
end if 
set H 



optimal 



: = n LP( Sensed U Optimal U jeAT ta (i) V] 



return (H t 



optimal ' sensed / 



Finally, we collect some straightforward facts about this algorithm; we omit the proof in the interest of brevity. 

Proposition VI.3 (Properties of the eight half-planes consensus algorithm) Consider a connected network Q sn of sensors 
that measure half-plane constraints and that implement the eight half-plane consensus algorithm. Assume the target does not 
move, that is, set v m . dx — 0. The following statements hold: 

(i) the candidate optimal constraints at each node contain the target at each instant of time; 

(ii) the candidate optimal constraints at each node monotonically improve over time; and 

(iii) additionally, if each node makes at most m measurements in finite time, then the candidate optimal constraints at each 
node converge in finite time to the globally optimal 8 half-plane constraints. 

Next, let us state some memory complexity bounds for the eight half-plane consensus algorithm. Again, we adopt the 
convention that a memory unit is the amount of memory required to store a constraint in H. Each node i of the network 
requires (8 + m + 8 card (A/in (*))) memory units in order to implement the algorithm. If we assume that number m of stored 
measurements is independent of n and that the indegree of each node is bounded irrespectively of n, then the algorithm memory 
complexity is in 0(1). Vice-versa, in worst-case graphs, the algorithm memory complexity is in 0(n). 

VII. Application to formation control for robotic networks 

In this section we apply constraints consensus ideas to formation control problems for networks of mobile robots. We focus 
on formations with the shapes of a point, a line, or a circle. (The problem of formation control to a point is usually referred to 
as the rendezvous or gathering problem.) We solve these formation control problems in a time-efficient manner via a distributed 
algorithm regulating the communication among robots and the motion of each robot. 



A. Model of robotic network 

We define a robotic network as follows. Each robot is equipped with a processor and robots exchange information via a 
communication graph. Therefore, the group of robots has the features of a synchronous network and can implement distributed 
algorithms as defined in Section III However, as compared with a synchronous network, a robotic network has two distinctions: 
(i) robots control their motion in space, and (ii) the communication graph among the robots depends upon the robots positions, 
rather than time. 

Specifically, the robotic network evolves according to the following discrete-time communication, computation and motion 
model. Each robot i £ {!,..., n} moves between rounds according to the first order discrete-time integrator p^(t + 1) = 
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P^(t) + vfi(t), where pM g R 2 and ||u^|| 2 < u mRX > 0. At each discrete time instant, robots at positions P n communicate 
according to the disk graph G<&sk(Pn) = ({1, • - • > n}, E^(P n )) defined as follows: an edge (i, j) € {1, . . . , n} 2 , i ^ j, belongs 
to E disk (P n ) if and only if ||p w - p b] \\ < r cmm for some r cmm > 0. 

A distributed algorithm for a robotic network consists of (1) a distributed algorithm for a synchronous network, that is, a 



processor state, a message alphabet, a message-generation and a state-transition function, as described in Section III (2) an 



additional function, called the control function, that determines the robot motion, with the following domain and co-domain: 

ctl : R 2 x W x A" -► B(0,u m . dx ). 

Additionally, we here allow the message generation and the state transition to depend upon not only the processor state but also 
the robot position. The state of the robotic network evolves as follows. First, at each communication round t, each processor 
i sends to its outgoing neighbors a message computed by applying the message-generation function to the current values of 
pW and W™. After a negligible period of time, the ith processor resets the value of its processor state iv™ by applying the 
state-transition function to the current values of pM and and to the messages received at time t. Finally, the position 
pM of the ith robot at time t is determined by applying the control function to the current value of pW and tt;M, and to the 
messages received at time t. 

In formal terms, if j/M (t) denotes the message vector received at time t by agent i (with y^ (t) being the message received 
from agent j), then the evolution is determined by 

yf ] (t)=m S g(pW(t-l),wW(t-l),i), 

ti;W(t)=stf( p W(t-l),«jW(t-l) ) tfW(t)), 

p W(t) =p[ i ](t-l) + ctl(p w (t-l),w [i] (f),y [<1 (f)), 

with the convention that msg(p bl (t - 1), w [j] (t - = null if E disk (p^(t ~ 1), . . . ,p [n] (t - 1)). 

B. Formation tasks and related optimization problems 

Numerous definitions of robot formation are considered in the multi-agent literature. Here we consider a somehow specific 
situation. Let Saints > Siines, and Secies be the set of points, lines and circles in the plane, respectively. We refer to these three 
sets as the shape sets. We aim to lead all robots in a network to a single element of one of the shape sets If S is a selected 
shape set, the formation task is achieved by the robotic network if there exists a time T € Z >0 such that for all t > T, all 
robots i G {1, .. . ,n} satisfy p^(t) € s for some element s 6 S. Specifically, the point-formation, or rendezvous task requires 
all connected robots to be at the same position, the line-formation task requires all connected robots to be on the same line, 
and the circle-formation task requires all connected robots to be on the same circle. 

We are interested in distributed algorithms that achieve such formation tasks optimally with respect to a suitable cost function. 
For the point-formation and line-formation tasks, we aim to minimize completion time, i.e., the time required by all robots to 
reach a common shape. For the circle-formation task, we aim to minimize the product between the time required to reach a 
common circle, and the diameter of the common circle. 

Remark VII. 1 (Circle formation) For the circle-formation problem we do not select the completion time as cost function, 
because of the following reasons. The centralized version of the minimum time circle-formation problem is equivalent to finding 
the minimum-width annulus containing the point-set. For arbitrary data sets, the minimum-width annulus has arbitrarily large 
minimum radius and bears similarities with the solution to the smallest stripe problem. For some configurations, all points 
are contained only in a small fraction of the minimum-width annulus; this is not the solution we envision when we consider 
moving robots in a circle formation. Therefore, we consider, instead, the smallest-area annulus. This cost function penalizes 
both the difference of the radiuses of the annulus (width of the annulus) and their sum. □ 

The key property of the minimum-time point-formation task, minimum-time line-formation task, and optimum circle-formation 
task is that their centralized versions are equivalent to finding the smallest ball, stripe and annulus, respectively, enclosing the 
n agents' initial positions. We state these equivalences in the following lemma without proof. 

Lemma VII.2 (Optimal shapes from geometric optimization) Given a set of distinct points {pi, . . . ,p„} c R 2 , consider 
the three optimization problems: 



mm max 

peSpoims ...,n} 



\Pj -P\\, 



min max distfo^.f), 

ees lin „ je{i,...,n} J 

min max dist(pj, c) • radius (c), 

eSScircles j'6{l,...,n} 
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where radius(c) denotes the radius of the circle c. These three optimization problems are equivalent to the smallest enclos- 
ing ball, the smallest enclosing stripe (for points in stripe -generic position), and the smallest enclosing annulus problem, 
respectively. 

Therefore, they are abstract optimization problems with combinatorial dimension 3, 5 and A, respectively. 

We conclude this section with some useful notation. We let target_set({pi, . . . ,p n }) denote the point, the line or the 
circle equidistant from the boundary of the smallest enclosing ball, stripe or annulus, respectively. 



C. Connectivity assumption, objective and strategy 

We assume that the robotic network is connected at initial time, i.e., that the graph Gd\sk(Pn(0)) is connected, and we aim 
to achieve the formation task while guaranteeing that the state-dependent communication graph remains connected during the 
evolution. The following connectivity maintenance strategy was originally proposed in |28] and a comprehensive discussion is 
in 1 41 1. The key idea is to restrict the allowable motion of each robot so as to preserve the existing edges in the communication 
graph. We present this idea in three steps. 

First, in a network with communication edges E'disk, if agents i and j are neighbors at time t e Z>o, then we require that 
their positions at time t + 1 belong to 



2 ■> 2 ^mm 

If all neighbors of agent i at time t are at locations QM(t) — , , . j ^ j 5 then the (convex) constraint set of agent i is 

*tpW(t),QW(i))= n B(^±^,\r m 

qe{qi,—,qi} 

Second, given p and q in M. 2 and a convex closed set QcK 2 with p E Q, we introduce the from-to-inside function, denoted 



by fti and illustrated in Figure 1 1 that computes the point in the closed segment [p, q] which is at the same time closest to q 
and inside Q. Formally, 

q, if q e Q, 

[p,q]ndQ, xfqiQ. 



fti(p, q, Q) 




Fig. 11. Illustration of the fti function; courtesy of the authors of [41]. 

Third and final, we assume that, independent of whatever control algorithm dictates the robots evolution, if p\^- get (t) denote 
the desired target positions of agent i at time t + 1, then we allow robot i to move towards that location only so far as the 
constraint set allows. This is encoded by: 

pW (t + 1) = fti(pH (t) , pM, et (i) , *(pW (t) , Q W (t))) . 



D. Move-to-consensus-shape strategy 

The minimum-time point-formation, minimum-time line-formation, and optimum circle formation tasks appear intractable 
in their general form due to the state-dependent communication constraints. To attack these problem, we search for an efficient 
strategy that converges to the optimal one when the upper bound on the robot speed u max goes to zero or, in other words, 
when information transmission tends to be infinitely faster than motion. To design such a strategy, we aim to reach consensus 
on the centralized solution to the problem, i.e., the optimal shape, and use the solution as a reference for the agents motion. 
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A simple sequential solution is as follows: first, the agents compute the optimal shape via constraints consensus, and then, 
when consensus is achieved, they move toward the closest point in the target shape (point, line or circle). An improved strategy 
allows concurrent execution of constraints consensus and motion: while the constraints consensus algorithm is running, each 
agent moves toward the estimated target position while maintaining connectivity of the communication graph. We first provide 
an informal description. 

Move-to-consensus-shape strategy: The processor state at each robot i consists of a set B^ of S candidate optimal 
constraints and a binary variable halt^ 6 {0, 1}. The set B^ is initialized to pW(0) and halt^ is initialized to 0. At 
each communication round, the processor performs the following tasks: (i) it transmits pM and to its neighbors 
and acquires its neighbors' candidate constraints and their current position; (ii) it runs an instance of the constraints 
consensus algorithm for the geometric optimization program of interest (smallest enclosing ball, stripe or annulus); 
if the constraints consensus halting condition is satisfied, it sets halt^ to 1; (iii) it computes a robot target position 
based on the current estimate of the optimal shape; (iv) it moves the robot towards the target position while respect 
input constraint and, if halv^ is still zero, enforcing connectivity with its current neighbors. 
Next we give a pseudo-code description. We let P„(0) = {p [1] (0), . . . ,p [n] (0)}- 



Problem data : A robotic network and a shape set 

Algorithm : Move-to-Consensus-Shape 
Message alphabet: A = R 2 U (R 2 ) s U {null} 

Processor state: B w C P n (0) with card(BW) = S 

halt H € {0,1} 
Physical state: p 1 ' 1 e I 2 

Initialization: P [l] := {p [l] (0), . . . ,P [l] {0)} 

halt 1 * 1 := 

function msg(pM, j\ 
l: return (pW,flW) 

function stf(pH , halt M , y) 

% executed by .node i, with y. ■ = (p^l Rb'I) •= msgfB^l i) 
l: 5 tmp := #0(0)} u BW £?( u;^ nW W) 

2: := Subex_LP(5 tmp ,pW) 

3: if B^ has not changed for 2n rounds, 

4: then halt 11 := 1; end if 

5: return (B^,halt [i] ) 

function ctl(B^,p^,y) 

% executed by node i, with yj = (jp&\ B^) := msg(B^,i) 

% map target_set defined in Section IV7/-BI 
1: Strgt := target_set(SW) 

2: ptrgt := arg min \\p [l] - p\\ 

PS Strgt 

3: if halt [l] = 0, then 

else X cnstr := B (pW , u max ) ; end if 
4: return fti(p^,p tIgU X cnstr ) -p® 

We refer the interested reader to [1| for numerical simulation results for the move-to-consensus-shape strategy. Finally, we 
state the correctness of the algorithm and omit the proof in the interest of brevity. 

Proposition VII.3 (Properties of the move-to-consensus-shape algorithm) On a robotic network with communication graph 
Gdi.sk and bounded control inputs u mm , the move-to-consensus-shape strategy achieves the desired formation control tasks. 

In the limit as u max — > + , the move-to-consensus-shape strategy solves the optimal formation control tasks, i.e., the 
minimum-time point-formation, minimum-time line-formation, and optimum circle formation tasks. 
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VIII. Conclusions 

We have introduced a novel class of distributed optimization problems and proposed a simple and intuitive algorithmic 
approach. Additionally, we have rigorously established the correctness of the proposed algorithms and presented a thorough 
Monte Carlo analysis to substantiate our conjecture that, for a broad variety of settings, the time complexity of our algorithms 
is linear. Finally, we have discussed in detail two modern applications: target localization in sensor networks and formation 
control in robotic networks. 

Promising avenues for further research include proving the linear-time-complexity conjecture, as well as applying our 
algorithms to (i) optimization problems for randomly switching graphs and gossip communication, (ii) distributed machine 
learning problems 11221 . Il2p . Il23ll . and (iii) convex quadratic programming [42|. Additionally, it is of interest to verify the 
performance of our proposed target localization and formation control algorithms in experimental setups. 

Appendix 

Abstract framework for the smallest enclosing stripe problem 

In this appendix we consider the smallest enclosing stripe problem in Example |H.2| for stripe-generic points and show that is 
satisfies the abstract optimization axioms. We begin with some notation. Let H be the set of constraints, i.e., the set of points 
in the plane for which we aim to compute the smallest enclosing stripe. With a slight abuse of notation, h € H denotes both 
a constraint and a point depending on the context. Let <fi(G), G C H, be the width of the smallest stripe enclosing the points 
in G. A constraint h is violated by the set G, i.e., <f>(G U {h}) > 4>(G), if the point h does not belong to the smallest stripe 
enclosing G. 

First, note that three points are necessary, but in general not sufficient, to uniquely identify the smallest enclosing stripe 
of a set H of three or more points. Indeed, the boundary of the smallest stripe enclosing H contains at least three points 
of H, e.g., point hi on one line and points {h2,hs} on the other line. So, the smallest enclosing stripe identifies the triplet 
{hi, h 2 , h 3 }. However, a simple geometric argument shows that the converse is not true, i.e., three points are not sufficient. 
Three points would be sufficient if the triplet was an ordered set and uniquely identified a stripe. To uniquely identify the 
smallest enclosing stripe for a set H of five or more points, one needs to add to the triplet {hi, h 2 , /13} another two points of 
H that belong to the correctly stripe among the three stripes determined by {hi, h 2 , h 3 }. In summary, this discussion shows 
that the combinatorial dimension is 5. 

Second, we prove that the two axioms of abstract optimization are satisfied. As usual, monotonicity is trivially satisfied, 
thus we need to prove locality. Suppose, by contradiction, that locality does not hold. Therefore there exist F, G and h, 
F <Z G <Z H and h <E H, with -00 < <f>(F) = 0(G) such that <f>(G U {h}) > 0(G) and <j>{F U {h}) = 0(F). Let S G and S F 
be the smallest stripes for G and F, so that (gi, 92,93) and (f\,f2,h) are the ordered triplets of points defining Sq and Sp, 
respectively. Since h is violated by G but not by F, this means that h belongs to Sf but not to Sg- Therefore, the stripes 
Sq and Sp must be different, i.e., (51,52,53) 7^ (fi,fi,h)- Because the points are in stripe-generic positions, we know that 
dist(<jri, I (<?2)<?3 )) 7^ dist(/i, l{fi, hj) and, therefore, (j>(G) > 0(F). This proves the contradiction. 
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